library(brms)
library(ggplot2)
library(lme4)
library(performance)
library(see)
library(sjmisc)
library(tidyverse)
options(
digits = 3
)
set.seed(170819)Analyses
Set-up
Packages
Data
d <- read_csv("data/data.csv")Rows: 960 Columns: 464
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (45): persistence, anonymity, topic, QUESTNNR_T1, MODE_T1, DE03_09_T1,...
dbl (407): case, n_OE, n_Words, mean_Words, number_messages, CASE_T1, REF_T...
lgl (8): SERIAL_T1, IV01_01_T1, K001_CP_T1, K001_T1, MAILSENT_T1, SERIAL_...
dttm (4): STARTED_T1, LASTDATA_T1, STARTED_T2, LASTDATA_T2
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.# same as above; but original file wording
# d <- read_csv("data/DataAggregated_T1T2_costsbenefits.csv")
# load("data/image.RData")
# recode to make as sum coding
d$anonymity_dev <- factor(d$anonymity)
contrasts(d$anonymity_dev) <- contr.sum(2)
d$persistence_dev <- factor(d$persistence)
contrasts(d$persistence_dev) <- contr.sum(2)
d <- d %>%
rename(
group = roles,
op_expr = n_OE
)Descriptives
Let’s inspect distribution of opinion expressions.
ggplot(d, aes(op_expr)) +
geom_histogram(binwidth = 1)
Looks like a zero-inflated poisson distribution. Confirms our preregistered approach to analyze data using zero-inflated Poisson approach.
nrow(d[d$op_expr == 0, ]) / nrow(d)[1] 0.214Overall, 21% of participants without any opinion expressions.
Let’s look at distribution of experimental groups.
d %>%
select(persistence, anonymity) %>%
table anonymity
persistence high low
high 240 240
low 240 240Distribution among groups perfect.
d %>%
select(topic) %>%
tabletopic
climate gender migration
320 320 320 Distribution of topics also perfect.
Let’s check if groups are nested in topics:
is_nested(d$topic, d$group)'f2' is nested within 'f1'[1] TRUEIndeed the case.
We first look at the experimental group’s descriptives
d %>%
group_by(persistence) %>%
summarize(op_expr_m = mean(op_expr))# A tibble: 2 × 2
persistence op_expr_m
<chr> <dbl>
1 high 9.26
2 low 9.29Looking at persistence, we see there’s virtually no difference among groups.
d %>%
group_by(anonymity) %>%
summarize(op_expr_m = mean(op_expr))# A tibble: 2 × 2
anonymity op_expr_m
<chr> <dbl>
1 high 9.03
2 low 9.51People who with less anonymity communicated more. But the difference isn’t particularly large.
d %>%
group_by(persistence, anonymity) %>%
summarize(op_expr_m = mean(op_expr))`summarise()` has grouped output by 'persistence'. You can override using the
`.groups` argument.# A tibble: 4 × 3
# Groups: persistence [2]
persistence anonymity op_expr_m
<chr> <chr> <dbl>
1 high high 9.27
2 high low 9.25
3 low high 8.79
4 low low 9.78Looking at both groups combined, we see that low anonymity and low persistence created highest participation. But differences among groups aren’t large.
d %>%
group_by(group) %>%
summarize(
anonymity = anonymity[1],
persistence = persistence[1],
topic = topic[1],
op_expr_m = mean(op_expr)
) %>%
rmarkdown::paged_table()Looking at the various individual groups, we do see some difference. Generally, this shows that individually, communication varied. But not so much systematically across experimental groups.
d %>%
group_by(topic) %>%
summarize(op_expr_m = mean(op_expr))# A tibble: 3 × 2
topic op_expr_m
<chr> <dbl>
1 climate 9.63
2 gender 9.96
3 migration 8.22Looking at topics specifically, we also see that there’s some variance.
Manipulation Check
Let’s see if respondents indeed perceived the experimental manipulations.
model_pers <- lm(
per_persistence ~ persistence_dev,
d
)
summary(model_pers)
Call:
lm(formula = per_persistence ~ persistence_dev, data = d)
Residuals:
Min 1Q Median 3Q Max
-3.194 -0.777 -0.110 0.806 3.223
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.9854 0.0378 79.1 <2e-16 ***
persistence_dev1 -1.2085 0.0378 -32.0 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 1 on 703 degrees of freedom
(255 observations deleted due to missingness)
Multiple R-squared: 0.593, Adjusted R-squared: 0.592
F-statistic: 1.02e+03 on 1 and 703 DF, p-value: <2e-16# report::report_text(model_pers)The experimental manipulation worked.
model_anon <- lm(
per_anonymity ~ anonymity_dev,
d
)
summary(model_anon)
Call:
lm(formula = per_anonymity ~ anonymity_dev, data = d)
Residuals:
Min 1Q Median 3Q Max
-3.425 -0.216 -0.216 0.575 3.784
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.8204 0.0348 81.0 <2e-16 ***
anonymity_dev1 -1.6042 0.0348 -46.1 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.924 on 703 degrees of freedom
(255 observations deleted due to missingness)
Multiple R-squared: 0.751, Adjusted R-squared: 0.751
F-statistic: 2.12e+03 on 1 and 703 DF, p-value: <2e-16# report::report_text(model_anon)The experimental manipulation worked.
Bayesian mixed effects modeling
We analyze the data using Bayesian modelling.
We use deviation/sum contrast coding (-.1, .1). Meaning, contrasts measure main effects of independent variables.
Fixed effects
We preregistered to analyze fixed effects. However, we did not explicate if we were to model zero inflation separately. Below hence two options.
fit_fe_1 <-
brm(
op_expr ~
1 + persistence_dev * anonymity_dev +
(1 | topic/group)
, data = d
, chains = 4
, cores = 4
, iter = 8000
, warmup = 2000
, family = zero_inflated_poisson()
, control = list(
adapt_delta = .95
, max_treedepth = 12
)
, save_pars = save_pars(all = TRUE)
)Running /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
using C compiler: ‘Apple clang version 15.0.0 (clang-1500.3.9.4)’
using SDK: ‘MacOSX14.4.sdk’
clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
In file included from <built-in>:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
#include <cmath>
^~~~~~~
1 error generated.
make: *** [foo.o] Error 1Warning: There were 767 divergent transitions after warmup. See
https://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
to find out why this is a problem and how to eliminate them.Warning: Examine the pairs() plot to diagnose sampling problemsWarning: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable.
Running the chains for more iterations may help. See
https://mc-stan.org/misc/warnings.html#tail-essOnly some convergence issues. Let’s inspect model.
plot(fit_fe_1, ask = FALSE)
Looks quite good!
Let’s look at results.
summary(fit_fe_1)Warning: There were 767 divergent transitions after warmup. Increasing
adapt_delta above 0.95 may help. See
http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup Family: zero_inflated_poisson
Links: mu = log; zi = identity
Formula: op_expr ~ 1 + persistence_dev * anonymity_dev + (1 | topic/group)
Data: d (Number of observations: 960)
Draws: 4 chains, each with iter = 8000; warmup = 2000; thin = 1;
total post-warmup draws = 24000
Multilevel Hyperparameters:
~topic (Number of levels: 3)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.26 0.28 0.01 1.05 1.00 1615 753
~topic:group (Number of levels: 48)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.40 0.05 0.32 0.51 1.00 2590 6421
Regression Coefficients:
Estimate Est.Error l-95% CI u-95% CI Rhat
Intercept 2.38 0.19 1.93 2.78 1.00
persistence_dev1 -0.00 0.06 -0.12 0.13 1.01
anonymity_dev1 0.00 0.06 -0.11 0.12 1.00
persistence_dev1:anonymity_dev1 0.03 0.06 -0.08 0.15 1.00
Bulk_ESS Tail_ESS
Intercept 1318 344
persistence_dev1 1207 369
anonymity_dev1 3229 6017
persistence_dev1:anonymity_dev1 3081 6134
Further Distributional Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
zi 0.21 0.01 0.19 0.24 1.00 6898 6452
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).No significant effect emerged.
Let’s inspect ICC
var_ratio_fe <- performance::variance_decomposition(
fit_fe_1
, by_group = TRUE)
var_ratio_fe# Random Effect Variances and ICC
Conditioned on: all random effects
## Variance Ratio (comparable to ICC)
Ratio: 0.43 CI 95%: [-0.14 0.73]
## Variances of Posterior Predicted Distribution
Conditioned on fixed effects: 29.50 CI 95%: [13.80 58.84]
Conditioned on rand. effects: 51.60 CI 95%: [46.86 56.86]
## Difference in Variances
Difference: 22.06 CI 95%: [-7.11 38.06]42.8177583 percent of variance in opinion expressions explained by both topics and groups.
Let’s look at exponentiated results for interpretation.
broom.mixed::tidy(fit_fe_1) |>
mutate(
estimate_exp = exp(estimate),
conf_low_exp = exp(conf.low),
conf_high_exp = exp(conf.high)
) %>%
mutate(
across(
c(
estimate
, conf.low
, conf.high
, estimate_exp
, conf_low_exp
, conf_high_exp
),
~ round(.x, 2)
)
) %>%
select(
-effect, -component, -group, -std.error
) %>%
rmarkdown::paged_table()Warning in tidy.brmsfit(fit_fe_1): some parameter names contain underscores:
term naming may be unreliable!Again, of course no significant effect. Main effects appear trivial. Interaction effect somewhat larger. Let’s visualize results to see what this exactly means.
p <- plot(
conditional_effects(
fit_fe_1
),
ask = FALSE,
plot = FALSE
)
p_anon <- p[["anonymity_dev"]]
p_int <- p[[3]]
p_pers <-
p[["persistence_dev"]] +
xlab("Persistence") +
ylab("Opinion expression") +
scale_x_discrete(
breaks = c(-0.5, .5),
labels = c("Low", "High")
) +
scale_y_continuous(
limits = c(5, 15)
)
p_anon <-
p[["anonymity_dev"]] +
xlab("Anonymity") +
ylab("Opinion expression") +
scale_x_discrete(
breaks = c(-0.5, .5),
labels = c("Low", "High")
) +
theme(
axis.title.y = element_blank()
) +
scale_y_continuous(
limits = c(5, 15)
)
p_int <-
p[[3]] +
xlab("Persistence") +
ylab("Opinion expression") +
scale_x_discrete(
breaks = c(-0.5, .5),
labels = c("Low", "High")
) +
scale_color_discrete(
labels = c("Low", "High")
) +
guides(
fill = "none",
color = guide_legend(
title = "Anonymity"
)
) +
theme(
axis.title.y = element_blank()
) +
scale_y_continuous(
limits = c(5, 15)
)
plot <- cowplot::plot_grid(
p_pers, p_anon, p_int,
labels = c('A', 'B', "C"),
nrow = 1,
rel_widths = c(2, 2, 3)
)
plot
ggsave("figures/results.png", plot, width = 8, height = 4)Shows that there are no main effects. There seems to be a (nonsignificant) interaction effect. In high persistence environment, anonymity is conducive to communication; in low it’s the opposite.
Let’s look at posteriors
p_1 <-
pp_check(fit_fe_1) +
labs(title = "Zero-inflated poisson")Using 10 posterior draws for ppc type 'dens_overlay' by default.p_1
The actual distribution cannot be precisely reproduced, but it’s also not too far off.
Random effects
We said we’d explore random effects. Following the “keep it maximal” principle, it’d actually make sense to consider these analyses the best ones, as they can model how the experimental conditions affect the outcomes differentially depending on topic.
Let’s now model random effects.
fit_re_1 <-
brm(
op_expr ~
1 + persistence_dev * anonymity_dev +
(1 + persistence_dev * anonymity_dev | topic) +
(1 | topic:group)
, data = d
, chains = 4
, cores = 4
, iter = 8000
, warmup = 2000
, family = zero_inflated_poisson()
, control = list(
adapt_delta = .95
, max_treedepth = 15
)
, save_pars = save_pars(all = TRUE)
)Running /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
using C compiler: ‘Apple clang version 15.0.0 (clang-1500.3.9.4)’
using SDK: ‘MacOSX14.4.sdk’
clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
In file included from <built-in>:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
#include <cmath>
^~~~~~~
1 error generated.
make: *** [foo.o] Error 1Warning: There were 1281 divergent transitions after warmup. See
https://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
to find out why this is a problem and how to eliminate them.Warning: Examine the pairs() plot to diagnose sampling problemsShows some problems, but not all that many.
Let’s inspect model.
plot(fit_re_1, ask = FALSE)
Traceplots look alright.
Let’s look at results.
summary(fit_re_1)Warning: There were 1281 divergent transitions after warmup. Increasing
adapt_delta above 0.95 may help. See
http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup Family: zero_inflated_poisson
Links: mu = log; zi = identity
Formula: op_expr ~ 1 + persistence_dev * anonymity_dev + (1 + persistence_dev * anonymity_dev | topic) + (1 | topic:group)
Data: d (Number of observations: 960)
Draws: 4 chains, each with iter = 8000; warmup = 2000; thin = 1;
total post-warmup draws = 24000
Multilevel Hyperparameters:
~topic (Number of levels: 3)
Estimate Est.Error
sd(Intercept) 0.36 0.46
sd(persistence_dev1) 0.25 0.35
sd(anonymity_dev1) 0.24 0.34
sd(persistence_dev1:anonymity_dev1) 0.41 0.46
cor(Intercept,persistence_dev1) -0.00 0.46
cor(Intercept,anonymity_dev1) 0.01 0.46
cor(persistence_dev1,anonymity_dev1) -0.01 0.47
cor(Intercept,persistence_dev1:anonymity_dev1) 0.08 0.46
cor(persistence_dev1,persistence_dev1:anonymity_dev1) 0.02 0.46
cor(anonymity_dev1,persistence_dev1:anonymity_dev1) -0.01 0.46
l-95% CI u-95% CI Rhat
sd(Intercept) 0.01 1.60 1.00
sd(persistence_dev1) 0.01 1.21 1.00
sd(anonymity_dev1) 0.00 1.20 1.00
sd(persistence_dev1:anonymity_dev1) 0.02 1.69 1.00
cor(Intercept,persistence_dev1) -0.83 0.83 1.00
cor(Intercept,anonymity_dev1) -0.83 0.83 1.00
cor(persistence_dev1,anonymity_dev1) -0.83 0.83 1.00
cor(Intercept,persistence_dev1:anonymity_dev1) -0.78 0.86 1.00
cor(persistence_dev1,persistence_dev1:anonymity_dev1) -0.83 0.84 1.00
cor(anonymity_dev1,persistence_dev1:anonymity_dev1) -0.83 0.83 1.00
Bulk_ESS Tail_ESS
sd(Intercept) 6620 11752
sd(persistence_dev1) 7764 10211
sd(anonymity_dev1) 7280 9370
sd(persistence_dev1:anonymity_dev1) 6024 5585
cor(Intercept,persistence_dev1) 21881 16269
cor(Intercept,anonymity_dev1) 22853 15065
cor(persistence_dev1,anonymity_dev1) 19146 13364
cor(Intercept,persistence_dev1:anonymity_dev1) 22045 15452
cor(persistence_dev1,persistence_dev1:anonymity_dev1) 19117 15472
cor(anonymity_dev1,persistence_dev1:anonymity_dev1) 15719 17072
~topic:group (Number of levels: 48)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.40 0.05 0.32 0.51 1.00 7776 12613
Regression Coefficients:
Estimate Est.Error l-95% CI u-95% CI Rhat
Intercept 2.38 0.27 1.79 2.95 1.00
persistence_dev1 -0.00 0.21 -0.41 0.44 1.00
anonymity_dev1 0.00 0.22 -0.44 0.44 1.00
persistence_dev1:anonymity_dev1 0.02 0.33 -0.67 0.71 1.00
Bulk_ESS Tail_ESS
Intercept 10960 9291
persistence_dev1 10694 8366
anonymity_dev1 9075 6054
persistence_dev1:anonymity_dev1 8583 5583
Further Distributional Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
zi 0.21 0.01 0.19 0.24 1.00 32167 15762
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).Again, virtually no main effect. Interaction effect larger, but also not even closely significant. Let’s see if the random effects model fits better
fit_fe_1 <- add_criterion(fit_fe_1, "loo", moment_match = TRUE)Warning: Some Pareto k diagnostic values are too high. See help('pareto-k-diagnostic') for details.Warning: Found 13 observations with a pareto_k > 0.7 in model 'fit_fe_1'. With
this many problematic observations, it may be more appropriate to use 'kfold'
with argument 'K = 10' to perform 10-fold cross-validation rather than LOO.fit_re_1 <- add_criterion(fit_re_1, "loo", moment_match = TRUE)Error in validate_ll(log_ratios) :
All input values must be finite or -Inf.
Error in eval(expr, envir, enclos) :
Exception: model6c7925fe16d3__namespace::write_array: Cor_1 is not positive definite. (in 'anon_model', line 188, column 2 to column 66)
Error in eval(expr, envir, enclos) :
Exception: model6c7925fe16d3__namespace::write_array: Cor_1 is not positive definite. (in 'anon_model', line 188, column 2 to column 66)Warning: Some Pareto k diagnostic values are too high. See help('pareto-k-diagnostic') for details.Warning: Found 17 observations with a pareto_k > 0.7 in model 'fit_re_1'. With
this many problematic observations, it may be more appropriate to use 'kfold'
with argument 'K = 10' to perform 10-fold cross-validation rather than LOO.comp_1 <- loo_compare(fit_fe_1, fit_re_1, criterion = "loo")
comp_1 elpd_diff se_diff
fit_fe_1 0.0 0.0
fit_re_1 -1.0 2.4 Although model comparisons showed that the model with random effects fitted better, the difference was not significant (Δ ELPD = -0.9779785, 95% CI [-5.7077202, 3.7517632]. Hence, for reasons of parsimony the model with fixed effects is preferred.
Hurdle
Let’s now estimate a fixed effects model with hurdles.
fit_hrdl_1 <-
brm(
bf(
op_expr ~
1 + persistence_dev * anonymity_dev +
(1 | topic) +
(1 | topic:group),
zi ~
1 + persistence_dev * anonymity_dev +
(1 | topic) +
(1 | topic:group)
)
, data = d
, chains = 4
, cores = 4
, iter = 8000
, warmup = 2000
, family = zero_inflated_poisson()
, control = list(
adapt_delta = .95
, max_treedepth = 15
)
)Running /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
using C compiler: ‘Apple clang version 15.0.0 (clang-1500.3.9.4)’
using SDK: ‘MacOSX14.4.sdk’
clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
In file included from <built-in>:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
#include <cmath>
^~~~~~~
1 error generated.
make: *** [foo.o] Error 1Warning: There were 601 divergent transitions after warmup. See
https://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
to find out why this is a problem and how to eliminate them.Warning: Examine the pairs() plot to diagnose sampling problemsEstimation works quite well.
Let’s inspect model.
plot(fit_hrdl_1, ask = FALSE)
Trace-plots look alright.
summary(fit_hrdl_1)Warning: There were 601 divergent transitions after warmup. Increasing
adapt_delta above 0.95 may help. See
http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup Family: zero_inflated_poisson
Links: mu = log; zi = logit
Formula: op_expr ~ 1 + persistence_dev * anonymity_dev + (1 | topic) + (1 | topic:group)
zi ~ 1 + persistence_dev * anonymity_dev + (1 | topic) + (1 | topic:group)
Data: d (Number of observations: 960)
Draws: 4 chains, each with iter = 8000; warmup = 2000; thin = 1;
total post-warmup draws = 24000
Multilevel Hyperparameters:
~topic (Number of levels: 3)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.27 0.32 0.01 1.17 1.00 2245 1183
sd(zi_Intercept) 0.28 0.44 0.01 1.50 1.00 6240 5389
~topic:group (Number of levels: 48)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.40 0.05 0.32 0.50 1.00 5272 10476
sd(zi_Intercept) 0.24 0.14 0.01 0.53 1.00 5361 8349
Regression Coefficients:
Estimate Est.Error l-95% CI u-95% CI Rhat
Intercept 2.39 0.20 1.96 2.83 1.00
zi_Intercept -1.30 0.28 -1.74 -0.67 1.00
persistence_dev1 -0.01 0.06 -0.12 0.11 1.00
anonymity_dev1 0.00 0.06 -0.12 0.12 1.00
persistence_dev1:anonymity_dev1 0.03 0.06 -0.09 0.15 1.00
zi_persistence_dev1 0.03 0.09 -0.15 0.20 1.00
zi_anonymity_dev1 0.01 0.09 -0.17 0.18 1.00
zi_persistence_dev1:anonymity_dev1 0.01 0.09 -0.16 0.18 1.00
Bulk_ESS Tail_ESS
Intercept 2994 1307
zi_Intercept 7245 4630
persistence_dev1 4984 9340
anonymity_dev1 5664 9023
persistence_dev1:anonymity_dev1 5982 9817
zi_persistence_dev1 19443 9630
zi_anonymity_dev1 22174 16988
zi_persistence_dev1:anonymity_dev1 10375 7268
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).Same results, no main effects, slightly larger but still nonsignificant interaction effect.
Let’s inspect coefficients individually.
coefficients(fit_hrdl_1)$topic
, , Intercept
Estimate Est.Error Q2.5 Q97.5
climate 2.441330 0.09198736 2.268439 2.630276
gender 2.416212 0.09056169 2.244328 2.603659
migration 2.306579 0.09761436 2.106201 2.485353
, , zi_Intercept
Estimate Est.Error Q2.5 Q97.5
climate -1.312747 0.1255348 -1.561256 -1.056106
gender -1.331712 0.1254764 -1.587955 -1.087517
migration -1.340756 0.1259943 -1.602138 -1.098537
, , persistence_dev1
Estimate Est.Error Q2.5 Q97.5
climate -0.005103491 0.05928069 -0.1215837 0.1109697
gender -0.005103491 0.05928069 -0.1215837 0.1109697
migration -0.005103491 0.05928069 -0.1215837 0.1109697
, , anonymity_dev1
Estimate Est.Error Q2.5 Q97.5
climate 0.0009798597 0.05871981 -0.1159714 0.1169367
gender 0.0009798597 0.05871981 -0.1159714 0.1169367
migration 0.0009798597 0.05871981 -0.1159714 0.1169367
, , persistence_dev1:anonymity_dev1
Estimate Est.Error Q2.5 Q97.5
climate 0.02947201 0.05854357 -0.08580015 0.1467483
gender 0.02947201 0.05854357 -0.08580015 0.1467483
migration 0.02947201 0.05854357 -0.08580015 0.1467483
, , zi_persistence_dev1
Estimate Est.Error Q2.5 Q97.5
climate 0.02959441 0.08905165 -0.148826 0.2034417
gender 0.02959441 0.08905165 -0.148826 0.2034417
migration 0.02959441 0.08905165 -0.148826 0.2034417
, , zi_anonymity_dev1
Estimate Est.Error Q2.5 Q97.5
climate 0.006613642 0.08968717 -0.1673214 0.1817552
gender 0.006613642 0.08968717 -0.1673214 0.1817552
migration 0.006613642 0.08968717 -0.1673214 0.1817552
, , zi_persistence_dev1:anonymity_dev1
Estimate Est.Error Q2.5 Q97.5
climate 0.008806348 0.08889244 -0.1644127 0.1844408
gender 0.008806348 0.08889244 -0.1644127 0.1844408
migration 0.008806348 0.08889244 -0.1644127 0.1844408
$`topic:group`
, , Intercept
Estimate Est.Error Q2.5 Q97.5
climate_Madrid1 1.759102 0.2462320 1.239356 2.240928
climate_Madrid2 2.457696 0.2344562 1.953004 2.931168
climate_Madrid3 2.272688 0.2356103 1.761028 2.755240
climate_Madrid4 2.792250 0.2333683 2.290435 3.268066
climate_Oslo1 1.984540 0.2399568 1.473012 2.458905
climate_Oslo2 2.675209 0.2333150 2.186113 3.147148
climate_Oslo3 2.351983 0.2385897 1.854177 2.831614
climate_Oslo4 2.763715 0.2329267 2.260392 3.241184
climate_Prag1 2.238111 0.2407547 1.745460 2.709367
climate_Prag2 2.285453 0.2405552 1.779591 2.751184
climate_Prag3 2.799090 0.2332121 2.307405 3.258438
climate_Prag4 2.562857 0.2378152 2.068329 3.028966
climate_Rio1 2.393871 0.2370286 1.879201 2.869942
climate_Rio2 2.526772 0.2375326 2.019074 2.988996
climate_Rio3 2.564091 0.2333346 2.056815 3.018157
climate_Rio4 2.342742 0.2342666 1.835705 2.821173
gender_Dubai1 1.908363 0.2450134 1.391224 2.399666
gender_Dubai2 2.309364 0.2340309 1.815465 2.780017
gender_Dubai3 2.126913 0.2390267 1.631461 2.595698
gender_Dubai4 2.483903 0.2324739 1.993163 2.956572
gender_London1 2.295351 0.2327935 1.816049 2.776759
gender_London2 2.530184 0.2319773 2.054167 3.003108
gender_London3 3.291772 0.2277683 2.821858 3.762740
gender_London4 1.927815 0.2363204 1.441430 2.422490
gender_Rom1 2.478452 0.2313663 1.988047 2.964336
gender_Rom2 2.968613 0.2261618 2.482351 3.449909
gender_Rom3 2.114065 0.2429499 1.594716 2.612889
gender_Rom4 1.870895 0.2389636 1.378278 2.365217
gender_Tokio1 1.874823 0.2383897 1.370670 2.357965
gender_Tokio2 2.406254 0.2307610 1.914750 2.877738
gender_Tokio3 3.084374 0.2281562 2.598814 3.559181
gender_Tokio4 2.796354 0.2287945 2.301066 3.281697
migration_Berlin1 2.236043 0.2427117 1.775399 2.772168
migration_Berlin2 1.975575 0.2425887 1.512872 2.530394
migration_Berlin3 2.421517 0.2411588 1.957223 2.961730
migration_Berlin4 3.316509 0.2327748 2.875410 3.844989
migration_Florenz1 1.816609 0.2453584 1.330089 2.371305
migration_Florenz2 2.454008 0.2411424 1.982611 3.000546
migration_Florenz3 2.180384 0.2435751 1.720796 2.729374
migration_Florenz4 2.301909 0.2408945 1.850372 2.853787
migration_Paris1 1.973939 0.2450204 1.516412 2.509617
migration_Paris2 2.061176 0.2442268 1.598793 2.592717
migration_Paris3 2.271019 0.2386424 1.820523 2.806410
migration_Paris4 2.413627 0.2392892 1.963469 2.945058
migration_Sydney1 2.833173 0.2368622 2.389574 3.357629
migration_Sydney2 2.217617 0.2445668 1.761643 2.750194
migration_Sydney3 2.428997 0.2436732 1.970346 2.969473
migration_Sydney4 2.524944 0.2397629 2.073708 3.062403
, , zi_Intercept
Estimate Est.Error Q2.5 Q97.5
climate_Madrid1 -1.370307 0.3721058 -2.106880 -0.61528568
climate_Madrid2 -1.421178 0.3827236 -2.191532 -0.66295948
climate_Madrid3 -1.314318 0.3680023 -2.006572 -0.55018676
climate_Madrid4 -1.314233 0.3615676 -1.988889 -0.55378347
climate_Oslo1 -1.268489 0.3648985 -1.928773 -0.49054187
climate_Oslo2 -1.271029 0.3629712 -1.919742 -0.50455122
climate_Oslo3 -1.322019 0.3658883 -2.020575 -0.55728205
climate_Oslo4 -1.116809 0.3840140 -1.724377 -0.24611549
climate_Prag1 -1.308407 0.3614901 -1.990049 -0.55732678
climate_Prag2 -1.205797 0.3640367 -1.827430 -0.39745096
climate_Prag3 -1.306936 0.3632909 -1.982392 -0.53620838
climate_Prag4 -1.307129 0.3618862 -1.977822 -0.56085374
climate_Rio1 -1.253380 0.3629223 -1.907779 -0.47604325
climate_Rio2 -1.149052 0.3774966 -1.772990 -0.29298313
climate_Rio3 -1.473682 0.4021293 -2.313559 -0.73098405
climate_Rio4 -1.307096 0.3653043 -1.987188 -0.54806953
gender_Dubai1 -1.098982 0.3848154 -1.695485 -0.21730384
gender_Dubai2 -1.355206 0.3690983 -2.070361 -0.61014911
gender_Dubai3 -1.356229 0.3746776 -2.082277 -0.58732582
gender_Dubai4 -1.356196 0.3705898 -2.089175 -0.60571297
gender_London1 -1.412310 0.3804841 -2.182613 -0.68963756
gender_London2 -1.303418 0.3665494 -1.993471 -0.54567104
gender_London3 -1.305189 0.3662344 -1.990580 -0.53238568
gender_London4 -1.361289 0.3683366 -2.084745 -0.60160100
gender_Rom1 -1.417191 0.3798392 -2.191277 -0.68170207
gender_Rom2 -1.207867 0.3660861 -1.834885 -0.38533330
gender_Rom3 -1.006564 0.4178990 -1.622727 -0.05748882
gender_Rom4 -1.259744 0.3649339 -1.908092 -0.47301603
gender_Tokio1 -1.422424 0.3831665 -2.210404 -0.68507654
gender_Tokio2 -1.316292 0.3642416 -2.005070 -0.55994924
gender_Tokio3 -1.211316 0.3647677 -1.831489 -0.40453563
gender_Tokio4 -1.368475 0.3721193 -2.108016 -0.61569715
migration_Berlin1 -1.362900 0.3711355 -2.088316 -0.60658403
migration_Berlin2 -1.365905 0.3711943 -2.094868 -0.63066407
migration_Berlin3 -1.105907 0.3877544 -1.705316 -0.20242532
migration_Berlin4 -1.412864 0.3808719 -2.187128 -0.67527560
migration_Florenz1 -1.373810 0.3775674 -2.153421 -0.61702415
migration_Florenz2 -1.262593 0.3630209 -1.923779 -0.48091321
migration_Florenz3 -1.312413 0.3633158 -1.984073 -0.55945309
migration_Florenz4 -1.316648 0.3661400 -2.010891 -0.55513195
migration_Paris1 -1.255209 0.3621782 -1.901247 -0.48386193
migration_Paris2 -1.147654 0.3765091 -1.753525 -0.29559249
migration_Paris3 -1.355056 0.3703160 -2.089092 -0.58375085
migration_Paris4 -1.304116 0.3643284 -1.981056 -0.54851894
migration_Sydney1 -1.410283 0.3816626 -2.182991 -0.66806776
migration_Sydney2 -1.145921 0.3780903 -1.764725 -0.28776382
migration_Sydney3 -1.250684 0.3656062 -1.900912 -0.44764680
migration_Sydney4 -1.409058 0.3825074 -2.186340 -0.64406722
, , persistence_dev1
Estimate Est.Error Q2.5 Q97.5
climate_Madrid1 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Madrid2 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Madrid3 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Madrid4 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Oslo1 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Oslo2 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Oslo3 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Oslo4 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Prag1 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Prag2 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Prag3 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Prag4 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Rio1 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Rio2 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Rio3 -0.005103491 0.05928069 -0.1215837 0.1109697
climate_Rio4 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Dubai1 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Dubai2 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Dubai3 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Dubai4 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_London1 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_London2 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_London3 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_London4 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Rom1 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Rom2 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Rom3 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Rom4 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Tokio1 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Tokio2 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Tokio3 -0.005103491 0.05928069 -0.1215837 0.1109697
gender_Tokio4 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Berlin1 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Berlin2 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Berlin3 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Berlin4 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Florenz1 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Florenz2 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Florenz3 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Florenz4 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Paris1 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Paris2 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Paris3 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Paris4 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Sydney1 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Sydney2 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Sydney3 -0.005103491 0.05928069 -0.1215837 0.1109697
migration_Sydney4 -0.005103491 0.05928069 -0.1215837 0.1109697
, , anonymity_dev1
Estimate Est.Error Q2.5 Q97.5
climate_Madrid1 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Madrid2 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Madrid3 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Madrid4 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Oslo1 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Oslo2 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Oslo3 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Oslo4 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Prag1 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Prag2 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Prag3 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Prag4 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Rio1 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Rio2 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Rio3 0.0009798597 0.05871981 -0.1159714 0.1169367
climate_Rio4 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Dubai1 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Dubai2 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Dubai3 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Dubai4 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_London1 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_London2 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_London3 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_London4 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Rom1 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Rom2 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Rom3 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Rom4 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Tokio1 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Tokio2 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Tokio3 0.0009798597 0.05871981 -0.1159714 0.1169367
gender_Tokio4 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Berlin1 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Berlin2 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Berlin3 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Berlin4 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Florenz1 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Florenz2 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Florenz3 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Florenz4 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Paris1 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Paris2 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Paris3 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Paris4 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Sydney1 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Sydney2 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Sydney3 0.0009798597 0.05871981 -0.1159714 0.1169367
migration_Sydney4 0.0009798597 0.05871981 -0.1159714 0.1169367
, , persistence_dev1:anonymity_dev1
Estimate Est.Error Q2.5 Q97.5
climate_Madrid1 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Madrid2 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Madrid3 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Madrid4 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Oslo1 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Oslo2 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Oslo3 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Oslo4 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Prag1 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Prag2 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Prag3 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Prag4 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Rio1 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Rio2 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Rio3 0.02947201 0.05854357 -0.08580015 0.1467483
climate_Rio4 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Dubai1 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Dubai2 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Dubai3 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Dubai4 0.02947201 0.05854357 -0.08580015 0.1467483
gender_London1 0.02947201 0.05854357 -0.08580015 0.1467483
gender_London2 0.02947201 0.05854357 -0.08580015 0.1467483
gender_London3 0.02947201 0.05854357 -0.08580015 0.1467483
gender_London4 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Rom1 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Rom2 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Rom3 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Rom4 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Tokio1 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Tokio2 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Tokio3 0.02947201 0.05854357 -0.08580015 0.1467483
gender_Tokio4 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Berlin1 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Berlin2 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Berlin3 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Berlin4 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Florenz1 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Florenz2 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Florenz3 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Florenz4 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Paris1 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Paris2 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Paris3 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Paris4 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Sydney1 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Sydney2 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Sydney3 0.02947201 0.05854357 -0.08580015 0.1467483
migration_Sydney4 0.02947201 0.05854357 -0.08580015 0.1467483
, , zi_persistence_dev1
Estimate Est.Error Q2.5 Q97.5
climate_Madrid1 0.02959441 0.08905165 -0.148826 0.2034417
climate_Madrid2 0.02959441 0.08905165 -0.148826 0.2034417
climate_Madrid3 0.02959441 0.08905165 -0.148826 0.2034417
climate_Madrid4 0.02959441 0.08905165 -0.148826 0.2034417
climate_Oslo1 0.02959441 0.08905165 -0.148826 0.2034417
climate_Oslo2 0.02959441 0.08905165 -0.148826 0.2034417
climate_Oslo3 0.02959441 0.08905165 -0.148826 0.2034417
climate_Oslo4 0.02959441 0.08905165 -0.148826 0.2034417
climate_Prag1 0.02959441 0.08905165 -0.148826 0.2034417
climate_Prag2 0.02959441 0.08905165 -0.148826 0.2034417
climate_Prag3 0.02959441 0.08905165 -0.148826 0.2034417
climate_Prag4 0.02959441 0.08905165 -0.148826 0.2034417
climate_Rio1 0.02959441 0.08905165 -0.148826 0.2034417
climate_Rio2 0.02959441 0.08905165 -0.148826 0.2034417
climate_Rio3 0.02959441 0.08905165 -0.148826 0.2034417
climate_Rio4 0.02959441 0.08905165 -0.148826 0.2034417
gender_Dubai1 0.02959441 0.08905165 -0.148826 0.2034417
gender_Dubai2 0.02959441 0.08905165 -0.148826 0.2034417
gender_Dubai3 0.02959441 0.08905165 -0.148826 0.2034417
gender_Dubai4 0.02959441 0.08905165 -0.148826 0.2034417
gender_London1 0.02959441 0.08905165 -0.148826 0.2034417
gender_London2 0.02959441 0.08905165 -0.148826 0.2034417
gender_London3 0.02959441 0.08905165 -0.148826 0.2034417
gender_London4 0.02959441 0.08905165 -0.148826 0.2034417
gender_Rom1 0.02959441 0.08905165 -0.148826 0.2034417
gender_Rom2 0.02959441 0.08905165 -0.148826 0.2034417
gender_Rom3 0.02959441 0.08905165 -0.148826 0.2034417
gender_Rom4 0.02959441 0.08905165 -0.148826 0.2034417
gender_Tokio1 0.02959441 0.08905165 -0.148826 0.2034417
gender_Tokio2 0.02959441 0.08905165 -0.148826 0.2034417
gender_Tokio3 0.02959441 0.08905165 -0.148826 0.2034417
gender_Tokio4 0.02959441 0.08905165 -0.148826 0.2034417
migration_Berlin1 0.02959441 0.08905165 -0.148826 0.2034417
migration_Berlin2 0.02959441 0.08905165 -0.148826 0.2034417
migration_Berlin3 0.02959441 0.08905165 -0.148826 0.2034417
migration_Berlin4 0.02959441 0.08905165 -0.148826 0.2034417
migration_Florenz1 0.02959441 0.08905165 -0.148826 0.2034417
migration_Florenz2 0.02959441 0.08905165 -0.148826 0.2034417
migration_Florenz3 0.02959441 0.08905165 -0.148826 0.2034417
migration_Florenz4 0.02959441 0.08905165 -0.148826 0.2034417
migration_Paris1 0.02959441 0.08905165 -0.148826 0.2034417
migration_Paris2 0.02959441 0.08905165 -0.148826 0.2034417
migration_Paris3 0.02959441 0.08905165 -0.148826 0.2034417
migration_Paris4 0.02959441 0.08905165 -0.148826 0.2034417
migration_Sydney1 0.02959441 0.08905165 -0.148826 0.2034417
migration_Sydney2 0.02959441 0.08905165 -0.148826 0.2034417
migration_Sydney3 0.02959441 0.08905165 -0.148826 0.2034417
migration_Sydney4 0.02959441 0.08905165 -0.148826 0.2034417
, , zi_anonymity_dev1
Estimate Est.Error Q2.5 Q97.5
climate_Madrid1 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Madrid2 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Madrid3 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Madrid4 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Oslo1 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Oslo2 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Oslo3 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Oslo4 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Prag1 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Prag2 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Prag3 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Prag4 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Rio1 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Rio2 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Rio3 0.006613642 0.08968717 -0.1673214 0.1817552
climate_Rio4 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Dubai1 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Dubai2 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Dubai3 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Dubai4 0.006613642 0.08968717 -0.1673214 0.1817552
gender_London1 0.006613642 0.08968717 -0.1673214 0.1817552
gender_London2 0.006613642 0.08968717 -0.1673214 0.1817552
gender_London3 0.006613642 0.08968717 -0.1673214 0.1817552
gender_London4 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Rom1 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Rom2 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Rom3 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Rom4 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Tokio1 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Tokio2 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Tokio3 0.006613642 0.08968717 -0.1673214 0.1817552
gender_Tokio4 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Berlin1 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Berlin2 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Berlin3 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Berlin4 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Florenz1 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Florenz2 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Florenz3 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Florenz4 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Paris1 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Paris2 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Paris3 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Paris4 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Sydney1 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Sydney2 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Sydney3 0.006613642 0.08968717 -0.1673214 0.1817552
migration_Sydney4 0.006613642 0.08968717 -0.1673214 0.1817552
, , zi_persistence_dev1:anonymity_dev1
Estimate Est.Error Q2.5 Q97.5
climate_Madrid1 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Madrid2 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Madrid3 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Madrid4 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Oslo1 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Oslo2 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Oslo3 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Oslo4 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Prag1 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Prag2 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Prag3 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Prag4 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Rio1 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Rio2 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Rio3 0.008806348 0.08889244 -0.1644127 0.1844408
climate_Rio4 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Dubai1 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Dubai2 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Dubai3 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Dubai4 0.008806348 0.08889244 -0.1644127 0.1844408
gender_London1 0.008806348 0.08889244 -0.1644127 0.1844408
gender_London2 0.008806348 0.08889244 -0.1644127 0.1844408
gender_London3 0.008806348 0.08889244 -0.1644127 0.1844408
gender_London4 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Rom1 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Rom2 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Rom3 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Rom4 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Tokio1 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Tokio2 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Tokio3 0.008806348 0.08889244 -0.1644127 0.1844408
gender_Tokio4 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Berlin1 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Berlin2 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Berlin3 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Berlin4 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Florenz1 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Florenz2 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Florenz3 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Florenz4 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Paris1 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Paris2 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Paris3 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Paris4 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Sydney1 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Sydney2 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Sydney3 0.008806348 0.08889244 -0.1644127 0.1844408
migration_Sydney4 0.008806348 0.08889244 -0.1644127 0.1844408Let’s look at exponentiated results for interpretation.
broom.mixed::tidy(fit_hrdl_1) |>
mutate(
estimate_exp = exp(estimate),
conf_low_exp = exp(conf.low),
conf_high_exp = exp(conf.high),
estimate_zi = exp(estimate) / (1 + exp(estimate)),
conf_low_zi = exp(conf.low) / (1 + exp(conf.low)),
conf_high_zi = exp(conf.high) / (1 + exp(conf.high))
) %>%
mutate(
across(
c(
estimate
, conf.low
, conf.high
, estimate_exp
, conf_low_exp
, conf_high_exp
, estimate_zi
, conf_low_zi
, conf_high_zi
),
~ round(.x, 2)
)
) %>%
select(
-effect, -group, -std.error # -component,
) %>%
rmarkdown::paged_table()Warning in tidy.brmsfit(fit_hrdl_1): some parameter names contain underscores:
term naming may be unreliable!Let’s visualize results.
plot(
conditional_effects(
fit_hrdl_1
),
ask = FALSE
)
Similar picture. Effects appear even smaller.
Let’s look at posteriors
p_4 <-
pp_check(fit_hrdl_1) +
labs(title = "Zero-inflated poisson")Using 10 posterior draws for ppc type 'dens_overlay' by default.p_4
Exploratory Analyses
Frequentist
Look at results from a frequentist perspective.
Fixed effects
Estimate nested model.
fit_fe_1_frq <-
lmer(
op_expr ~
1 +
(1 | topic/group) +
persistence_dev * anonymity_dev
, data = d
)
summary(fit_fe_1_frq)Linear mixed model fit by REML ['lmerMod']
Formula: op_expr ~ 1 + (1 | topic/group) + persistence_dev * anonymity_dev
Data: d
REML criterion at convergence: 7606.1
Scaled residuals:
Min 1Q Median 3Q Max
-1.4108 -0.5454 -0.2147 0.2659 13.5411
Random effects:
Groups Name Variance Std.Dev.
group:topic (Intercept) 10.94 3.307
topic (Intercept) 0.00 0.000
Residual 155.90 12.486
Number of obs: 960, groups: group:topic, 48; topic, 3
Fixed effects:
Estimate Std. Error t value
(Intercept) 9.2729 0.6247 14.844
persistence_dev1 -0.0125 0.6247 -0.020
anonymity_dev1 -0.2417 0.6247 -0.387
persistence_dev1:anonymity_dev1 0.2521 0.6247 0.404
Correlation of Fixed Effects:
(Intr) prss_1 anny_1
prsstnc_dv1 0.000
annymty_dv1 0.000 0.000
prsstn_1:_1 0.000 0.000 0.000
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')Quite weird that topic doesn’t get any variance at all. Perhaps due to small cluster size? With Bayesian estimation, it worked alright.
Estimate without nesting.
fit_fe_2_frq <-
lmer(
op_expr ~
1 +
(1 | topic) +
persistence_dev * anonymity_dev
, data = d
)
summary(fit_fe_2_frq)Linear mixed model fit by REML ['lmerMod']
Formula: op_expr ~ 1 + (1 | topic) + persistence_dev * anonymity_dev
Data: d
REML criterion at convergence: 7627
Scaled residuals:
Min 1Q Median 3Q Max
-0.7806 -0.6106 -0.2378 0.2292 13.7610
Random effects:
Groups Name Variance Std.Dev.
topic (Intercept) 0.3331 0.5771
Residual 165.7407 12.8740
Number of obs: 960, groups: topic, 3
Fixed effects:
Estimate Std. Error t value
(Intercept) 9.2729 0.5326 17.410
persistence_dev1 -0.0125 0.4155 -0.030
anonymity_dev1 -0.2417 0.4155 -0.582
persistence_dev1:anonymity_dev1 0.2521 0.4155 0.607
Correlation of Fixed Effects:
(Intr) prss_1 anny_1
prsstnc_dv1 0.000
annymty_dv1 0.000 0.000
prsstn_1:_1 0.000 0.000 0.000 Estimate without hierarchical structure.
fit_fe_3_frq <-
lm(
op_expr ~
1 +
persistence_dev * anonymity_dev + topic
, data = d
)
summary(fit_fe_3_frq)
Call:
lm(formula = op_expr ~ 1 + persistence_dev * anonymity_dev +
topic, data = d)
Residuals:
Min 1Q Median 3Q Max
-10.469 -7.744 -3.150 3.040 177.798
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 9.6312 0.7197 13.383 <2e-16 ***
persistence_dev1 -0.0125 0.4155 -0.030 0.976
anonymity_dev1 -0.2417 0.4155 -0.582 0.561
topicgender 0.3312 1.0178 0.325 0.745
topicmigration -1.4062 1.0178 -1.382 0.167
persistence_dev1:anonymity_dev1 0.2521 0.4155 0.607 0.544
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 12.87 on 954 degrees of freedom
Multiple R-squared: 0.004169, Adjusted R-squared: -0.001051
F-statistic: 0.7987 on 5 and 954 DF, p-value: 0.5507Also here, no significant effects.
Random effects
Now do the same with random effects.
fit_re_1_frq <-
lmer(
op_expr ~
1 +
persistence_dev * anonymity_dev +
(1 + persistence_dev * anonymity_dev | topic) +
(1 | topic:group),
, data = d
)
summary(fit_re_1_frq)Linear mixed model fit by REML ['lmerMod']
Formula:
op_expr ~ 1 + persistence_dev * anonymity_dev + (1 + persistence_dev *
anonymity_dev | topic) + (1 | topic:group)
Data: d
REML criterion at convergence: 7603.3
Scaled residuals:
Min 1Q Median 3Q Max
-1.4310 -0.5485 -0.2174 0.2635 13.5713
Random effects:
Groups Name Variance Std.Dev. Corr
topic:group (Intercept) 8.751e+00 2.95816
topic (Intercept) 5.858e-01 0.76535
persistence_dev1 1.752e-02 0.13238 -1.00
anonymity_dev1 8.313e-03 0.09117 -1.00 1.00
persistence_dev1:anonymity_dev1 2.396e+00 1.54785 1.00 -1.00
Residual 1.559e+02 12.48581
-1.00
Number of obs: 960, groups: topic:group, 48; topic, 3
Fixed effects:
Estimate Std. Error t value
(Intercept) 9.2729 0.7348 12.619
persistence_dev1 -0.0125 0.5921 -0.021
anonymity_dev1 -0.2417 0.5895 -0.410
persistence_dev1:anonymity_dev1 0.2521 1.0693 0.236
Correlation of Fixed Effects:
(Intr) prss_1 anny_1
prsstnc_dv1 -0.078
annymty_dv1 -0.054 0.012
prsstn_1:_1 0.503 -0.108 -0.075
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')fit_re_2_frq <-
lmer(
op_expr ~
1 +
persistence_dev * anonymity_dev +
(1 + persistence_dev * anonymity_dev | topic),
, data = d
)
summary(fit_re_2_frq)Linear mixed model fit by REML ['lmerMod']
Formula:
op_expr ~ 1 + persistence_dev * anonymity_dev + (1 + persistence_dev *
anonymity_dev | topic)
Data: d
REML criterion at convergence: 7617.8
Scaled residuals:
Min 1Q Median 3Q Max
-0.9574 -0.5665 -0.2408 0.2263 13.7512
Random effects:
Groups Name Variance Std.Dev. Corr
topic (Intercept) 6.878e-01 0.82936
persistence_dev1 2.061e-02 0.14355 -1.00
anonymity_dev1 9.755e-03 0.09876 -1.00 1.00
persistence_dev1:anonymity_dev1 2.814e+00 1.67745 1.00 -1.00 -1.00
Residual 1.636e+02 12.79063
Number of obs: 960, groups: topic, 3
Fixed effects:
Estimate Std. Error t value
(Intercept) 9.2729 0.6322 14.667
persistence_dev1 -0.0125 0.4210 -0.030
anonymity_dev1 -0.2417 0.4167 -0.580
persistence_dev1:anonymity_dev1 0.2521 1.0528 0.239
Correlation of Fixed Effects:
(Intr) prss_1 anny_1
prsstnc_dv1 -0.149
annymty_dv1 -0.104 0.027
prsstn_1:_1 0.697 -0.181 -0.126
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')Let’s inspect coefficients.
coefficients(fit_re_2_frq)$topic
(Intercept) persistence_dev1 anonymity_dev1
climate 9.331303 -0.02260573 -0.2486197
gender 10.017024 -0.14129197 -0.3302795
migration 8.470422 0.12639770 -0.1461008
persistence_dev1:anonymity_dev1
climate 0.3701762
gender 1.7571129
migration -1.3710392
attr(,"class")
[1] "coef.mer"Gender
Let’s see if effects differ across genders
fit_fe_gen_1 <-
brm(
op_expr ~
1 + persistence_dev * anonymity_dev +
(1 | topic/group)
, data = d
, chains = 4
, cores = 4
, iter = 8000
, warmup = 2000
, family = zero_inflated_poisson()
, control = list(
adapt_delta = .95
, max_treedepth = 12
)
)Running /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
using C compiler: ‘Apple clang version 15.0.0 (clang-1500.3.9.4)’
using SDK: ‘MacOSX14.4.sdk’
clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
In file included from <built-in>:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
#include <cmath>
^~~~~~~
1 error generated.
make: *** [foo.o] Error 1Warning: There were 315 divergent transitions after warmup. See
https://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
to find out why this is a problem and how to eliminate them.Warning: Examine the pairs() plot to diagnose sampling problemsShows some problems, but not all that many.
Let’s inspect model.
plot(fit_fe_gen_1, ask = FALSE)
Traceplots look alright.
Let’s look at results.
summary(fit_fe_gen_1)Warning: There were 315 divergent transitions after warmup. Increasing
adapt_delta above 0.95 may help. See
http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup Family: zero_inflated_poisson
Links: mu = log; zi = identity
Formula: op_expr ~ 1 + persistence_dev * anonymity_dev + (1 | topic/group)
Data: d (Number of observations: 960)
Draws: 4 chains, each with iter = 8000; warmup = 2000; thin = 1;
total post-warmup draws = 24000
Multilevel Hyperparameters:
~topic (Number of levels: 3)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.27 0.33 0.01 1.21 1.00 2388 3643
~topic:group (Number of levels: 48)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.40 0.05 0.32 0.50 1.00 3576 5264
Regression Coefficients:
Estimate Est.Error l-95% CI u-95% CI Rhat
Intercept 2.39 0.20 1.95 2.82 1.00
persistence_dev1 -0.00 0.06 -0.12 0.11 1.00
anonymity_dev1 0.00 0.06 -0.12 0.12 1.00
persistence_dev1:anonymity_dev1 0.03 0.06 -0.09 0.15 1.00
Bulk_ESS Tail_ESS
Intercept 4302 2993
persistence_dev1 4061 6358
anonymity_dev1 4270 7001
persistence_dev1:anonymity_dev1 3856 5126
Further Distributional Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
zi 0.21 0.01 0.19 0.24 1.00 16100 13653
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).Benefits
Let’s see if benefits differ across experimental groups.
We first look at the experimental group’s descriptives
d %>%
group_by(persistence) %>%
summarize(benefits_m = mean(benefits, na.rm = TRUE))# A tibble: 2 × 2
persistence benefits_m
<chr> <dbl>
1 high 3.12
2 low 3.23Looking at persistence, we see people with lower persistence reporting slightly higher benefits.
d %>%
group_by(anonymity) %>%
summarize(benefits_m = mean(benefits, na.rm = TRUE))# A tibble: 2 × 2
anonymity benefits_m
<chr> <dbl>
1 high 3.14
2 low 3.20Looking at anonymity, we see people with low anonymity reporting marginally higher benefits.
d %>%
group_by(persistence, anonymity) %>%
summarize(benefits_m = mean(benefits, na.rm = T))`summarise()` has grouped output by 'persistence'. You can override using the
`.groups` argument.# A tibble: 4 × 3
# Groups: persistence [2]
persistence anonymity benefits_m
<chr> <chr> <dbl>
1 high high 3.07
2 high low 3.18
3 low high 3.22
4 low low 3.23Looking at both groups combined, we see that low anonymity and low persistence yielded highest benefits.
fit_fe_ben_1 <-
brm(
benefits ~
1 + persistence_dev * anonymity_dev +
(1 | topic/group)
, data = d
, chains = 4
, cores = 4
, iter = 8000
, warmup = 2000
, control = list(
adapt_delta = .95
, max_treedepth = 12
)
)Warning: Rows containing NAs were excluded from the model.Running /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
using C compiler: ‘Apple clang version 15.0.0 (clang-1500.3.9.4)’
using SDK: ‘MacOSX14.4.sdk’
clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
In file included from <built-in>:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
#include <cmath>
^~~~~~~
1 error generated.
make: *** [foo.o] Error 1Warning: There were 146 divergent transitions after warmup. See
https://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
to find out why this is a problem and how to eliminate them.Warning: Examine the pairs() plot to diagnose sampling problemsShows some problems, but not all that many.
Let’s inspect model.
plot(fit_fe_ben_1, ask = FALSE)
Traceplots look alright.
Let’s look at results.
summary(fit_fe_ben_1)Warning: There were 146 divergent transitions after warmup. Increasing
adapt_delta above 0.95 may help. See
http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup Family: gaussian
Links: mu = identity; sigma = identity
Formula: benefits ~ 1 + persistence_dev * anonymity_dev + (1 | topic/group)
Data: d (Number of observations: 705)
Draws: 4 chains, each with iter = 8000; warmup = 2000; thin = 1;
total post-warmup draws = 24000
Multilevel Hyperparameters:
~topic (Number of levels: 3)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.14 0.22 0.00 0.78 1.00 3274 2146
~topic:group (Number of levels: 48)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.07 0.05 0.00 0.17 1.00 6912 9608
Regression Coefficients:
Estimate Est.Error l-95% CI u-95% CI Rhat
Intercept 3.18 0.11 2.93 3.44 1.00
persistence_dev1 -0.05 0.03 -0.11 0.01 1.00
anonymity_dev1 -0.03 0.03 -0.09 0.03 1.00
persistence_dev1:anonymity_dev1 -0.02 0.03 -0.08 0.04 1.00
Bulk_ESS Tail_ESS
Intercept 5292 2048
persistence_dev1 25336 17342
anonymity_dev1 24256 16751
persistence_dev1:anonymity_dev1 24260 15633
Further Distributional Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sigma 0.74 0.02 0.70 0.78 1.00 29448 16416
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).No significant effect. But note that effect of persistence on perceived benefits only marginally not significant.
Risks
Let’s see if perceived differed across experimental groups.
We first look at the experimental group’s descriptives
d %>%
group_by(persistence) %>%
summarize(costs = mean(costs, na.rm = TRUE))# A tibble: 2 × 2
persistence costs
<chr> <dbl>
1 high 1.99
2 low 1.99Looking at persistence, we see both groups report equal costs.
d %>%
group_by(anonymity) %>%
summarize(costs = mean(costs, na.rm = TRUE))# A tibble: 2 × 2
anonymity costs
<chr> <dbl>
1 high 1.89
2 low 2.09Looking at anonymity, we see people with low anonymity report slightly higher costs.
d %>%
group_by(persistence, anonymity) %>%
summarize(costs = mean(costs, na.rm = TRUE))`summarise()` has grouped output by 'persistence'. You can override using the
`.groups` argument.# A tibble: 4 × 3
# Groups: persistence [2]
persistence anonymity costs
<chr> <chr> <dbl>
1 high high 1.90
2 high low 2.07
3 low high 1.87
4 low low 2.11Looking at both groups combined, we see that highest costs were reported by group with low anonymity and low persistence.
fit_fe_ris_1 <-
brm(
costs ~
1 + persistence_dev * anonymity_dev +
(1 | topic/group)
, data = d
, chains = 4
, cores = 4
, iter = 20000
, warmup = 2000
, control = list(
adapt_delta = .95
, max_treedepth = 12
)
)Warning: Rows containing NAs were excluded from the model.Running /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
using C compiler: ‘Apple clang version 15.0.0 (clang-1500.3.9.4)’
using SDK: ‘MacOSX14.4.sdk’
clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
In file included from <built-in>:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
#include <cmath>
^~~~~~~
1 error generated.
make: *** [foo.o] Error 1Warning: There were 495 divergent transitions after warmup. See
https://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
to find out why this is a problem and how to eliminate them.Warning: Examine the pairs() plot to diagnose sampling problemsShows some problems, but not all that many.
Let’s inspect model.
plot(fit_fe_ris_1, ask = FALSE)
Traceplots look alright.
Let’s look at results.
summary(fit_fe_ris_1)Warning: There were 495 divergent transitions after warmup. Increasing
adapt_delta above 0.95 may help. See
http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup Family: gaussian
Links: mu = identity; sigma = identity
Formula: costs ~ 1 + persistence_dev * anonymity_dev + (1 | topic/group)
Data: d (Number of observations: 705)
Draws: 4 chains, each with iter = 20000; warmup = 2000; thin = 1;
total post-warmup draws = 72000
Multilevel Hyperparameters:
~topic (Number of levels: 3)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.14 0.21 0.00 0.74 1.00 11201 5353
~topic:group (Number of levels: 48)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.07 0.05 0.00 0.17 1.00 19946 24053
Regression Coefficients:
Estimate Est.Error l-95% CI u-95% CI Rhat
Intercept 1.99 0.11 1.74 2.23 1.00
persistence_dev1 -0.00 0.03 -0.07 0.06 1.00
anonymity_dev1 -0.10 0.03 -0.17 -0.03 1.00
persistence_dev1:anonymity_dev1 0.02 0.03 -0.05 0.08 1.00
Bulk_ESS Tail_ESS
Intercept 10504 5540
persistence_dev1 62661 47237
anonymity_dev1 65846 48554
persistence_dev1:anonymity_dev1 58363 30821
Further Distributional Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sigma 0.85 0.02 0.81 0.90 1.00 75075 46472
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).We find that anonymity does reduce costs.
Mediation
Let’s see if perceived benefits and risks were associated with increased opinion expressions.
fit_fe_med <-
brm(
op_expr ~
1 + persistence_dev * anonymity_dev + benefits + costs +
(1 | topic/group)
, data = d
, chains = 4
, cores = 4
, iter = 8000
, warmup = 2000
, family = zero_inflated_poisson()
, control = list(
adapt_delta = .95
, max_treedepth = 12
)
)Warning: Rows containing NAs were excluded from the model.Compiling Stan program...Trying to compile a simple C fileRunning /Library/Frameworks/R.framework/Resources/bin/R CMD SHLIB foo.c
using C compiler: ‘Apple clang version 15.0.0 (clang-1500.3.9.4)’
using SDK: ‘MacOSX14.4.sdk’
clang -arch arm64 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/Rcpp/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/unsupported" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/BH/include" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/src/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppParallel/include/" -I"/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/rstan/include" -DEIGEN_NO_DEBUG -DBOOST_DISABLE_ASSERTS -DBOOST_PENDING_INTEGER_LOG2_HPP -DSTAN_THREADS -DUSE_STANC3 -DSTRICT_R_HEADERS -DBOOST_PHOENIX_NO_VARIADIC_EXPRESSION -D_HAS_AUTO_PTR_ETC=0 -include '/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp' -D_REENTRANT -DRCPP_PARALLEL_USE_TBB=1 -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c foo.c -o foo.o
In file included from <built-in>:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/StanHeaders/include/stan/math/prim/fun/Eigen.hpp:22:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Dense:1:
In file included from /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/Core:19:
/Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/library/RcppEigen/include/Eigen/src/Core/util/Macros.h:679:10: fatal error: 'cmath' file not found
#include <cmath>
^~~~~~~
1 error generated.
make: *** [foo.o] Error 1Start samplingWarning: There were 351 divergent transitions after warmup. See
https://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
to find out why this is a problem and how to eliminate them.Warning: Examine the pairs() plot to diagnose sampling problemsLet’s look at results.
summary(fit_fe_med)Warning: There were 351 divergent transitions after warmup. Increasing
adapt_delta above 0.95 may help. See
http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup Family: zero_inflated_poisson
Links: mu = log; zi = identity
Formula: op_expr ~ 1 + persistence_dev * anonymity_dev + benefits + costs + (1 | topic/group)
Data: d (Number of observations: 705)
Draws: 4 chains, each with iter = 8000; warmup = 2000; thin = 1;
total post-warmup draws = 24000
Multilevel Hyperparameters:
~topic (Number of levels: 3)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.21 0.25 0.01 0.96 1.00 3142 7447
~topic:group (Number of levels: 48)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.41 0.05 0.33 0.51 1.00 4402 7783
Regression Coefficients:
Estimate Est.Error l-95% CI u-95% CI Rhat
Intercept 2.31 0.16 1.94 2.66 1.00
persistence_dev1 -0.02 0.06 -0.14 0.10 1.00
anonymity_dev1 0.01 0.06 -0.11 0.13 1.00
benefits 0.11 0.02 0.08 0.14 1.00
costs -0.12 0.01 -0.14 -0.09 1.00
persistence_dev1:anonymity_dev1 0.03 0.06 -0.09 0.15 1.00
Bulk_ESS Tail_ESS
Intercept 6722 5749
persistence_dev1 3874 6750
anonymity_dev1 3995 6829
benefits 16490 11997
costs 16979 15573
persistence_dev1:anonymity_dev1 3881 7316
Further Distributional Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
zi 0.08 0.01 0.06 0.10 1.00 17823 13042
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).We find that increased perceived costs are associated with decreased opinion expressions. Increased benefits are associated with increased opinion expressions. Let’s check if overall effect is significant.
anon_risk_a_b <- fixef(fit_fe_ris_1)["anonymity_dev1", "Estimate"]
anon_risk_a_se <- fixef(fit_fe_ris_1)["anonymity_dev1", "Est.Error"]
anon_risk_a_dis <- rnorm(10000, anon_risk_a_b, anon_risk_a_se)
anon_risk_b_b <- fixef(fit_fe_med)["benefits", "Estimate"]
anon_risk_b_se <- fixef(fit_fe_med)["benefits", "Est.Error"]
anon_risk_b_dis <- rnorm(10000, anon_risk_b_b, anon_risk_b_se)
anon_risk_ab_dis <- anon_risk_a_dis * anon_risk_b_dis
anon_risk_ab_m <- median(anon_risk_ab_dis)
anon_risk_ab_ll <- quantile(anon_risk_ab_dis, .025)
anon_risk_ab_ul <- quantile(anon_risk_ab_dis, .975)The effect is significant (b = -0.0107919, 95% MC CI [-0.0198564, -0.0034263]).
Save
save.image("data/image.RData")